<!DOCTYPE HTML>
<html>
<head>
<meta charset="UTF-8">
<title>l1 trend filtering</title>
<link rel="canonical" href="http://cvxr.com/cvx/examples/time_series_analysis/html/l1_trend_filter_snp500.html">
<link rel="stylesheet" href="../../examples.css" type="text/css">
</head>
<body>
<div id="header">
<h1>l1 trend filtering</h1>
Jump to:&nbsp;&nbsp;&nbsp;&nbsp;
<a href="#source">Source code</a>&nbsp;&nbsp;&nbsp;&nbsp;
<a href="#output">Text output</a>
&nbsp;&nbsp;&nbsp;&nbsp;
<a href="#plots">Plots</a>
&nbsp;&nbsp;&nbsp;&nbsp;<a href="../../index.html">Library index</a>
</div>
<div id="content">
<a id="source"></a>
<pre class="codeinput">
<span class="comment">% Written for CVX by Kwangmoo Koh - 12/10/07</span>
<span class="comment">%</span>
<span class="comment">% The problem of estimating underlying trends in time series data arises in</span>
<span class="comment">% a variety of disciplines. The l1 trend filtering method produces trend</span>
<span class="comment">% estimates x that are piecewise linear from the time series y.</span>
<span class="comment">%</span>
<span class="comment">% The l1 trend estimation problem can be formulated as</span>
<span class="comment">%</span>
<span class="comment">%    minimize    (1/2)*||y-x||^2+lambda*||Dx||_1,</span>
<span class="comment">%</span>
<span class="comment">% with variable x , and problem data y and lambda, with lambda &gt;0.</span>
<span class="comment">% D is the second difference matrix, with rows [0... -1 2 -1 ...0]</span>
<span class="comment">%</span>
<span class="comment">% CVX is not optimized for the l1 trend filtering problem.</span>
<span class="comment">% For large problems, use l1_tf (www.stanford.edu/~boyd/l1_tf/).</span>

<span class="comment">% load time series data</span>
y = csvread(<span class="string">'snp500.txt'</span>); <span class="comment">% log price of snp500</span>
n = length(y);

<span class="comment">% form second difference matrix</span>
e = ones(n,1);
D = spdiags([e -2*e e], 0:2, n-2, n);

<span class="comment">% set regularization parameter</span>
lambda = 50;

<span class="comment">% solve l1 trend filtering problem</span>
cvx_begin
    variable <span class="string">x(n)</span>
    minimize( 0.5*sum_square(y-x)+lambda*norm(D*x,1) )
cvx_end

<span class="comment">% plot estimated trend with original signal</span>
figure(1);
plot(1:n,y,<span class="string">'k:'</span>,<span class="string">'LineWidth'</span>,1.0); hold <span class="string">on</span>;
plot(1:n,x,<span class="string">'b-'</span>,<span class="string">'LineWidth'</span>,2.0); hold <span class="string">off</span>;
xlabel(<span class="string">'date'</span>); ylabel(<span class="string">'log price'</span>);
</pre>
<a id="output"></a>
<pre class="codeoutput">
 
Calling sedumi: 5998 variables, 1999 equality constraints
------------------------------------------------------------
SeDuMi 1.21 by AdvOL, 2005-2008 and Jos F. Sturm, 1998-2003.
Alg = 2: xz-corrector, Adaptive Step-Differentiation, theta = 0.250, beta = 0.500
eqs m = 1999, order n = 3999, dim = 5999, blocks = 2
nnz(A) = 9991 + 1, nnz(ADA) = 9985, nnz(L) = 5992
Handling 1 + 1 dense columns.
 it :     b*y       gap    delta  rate   t/tP*  t/tD*   feas cg cg  prec
  0 :            3.92E+05 0.000
  1 :  -1.27E+01 2.07E+03 0.000 0.0053 0.9990 0.9990   1.04  1  1  1.0E+00
  2 :  -2.68E+00 4.81E+02 0.000 0.2326 0.9000 0.9000   1.00  1  1  9.9E-01
  3 :   7.05E-02 1.15E+02 0.000 0.2400 0.9000 0.9000   1.00  1  1  9.8E-01
  4 :   6.04E-01 5.10E+01 0.000 0.4419 0.9000 0.9000   1.00  1  1  9.7E-01
  5 :   7.94E-01 2.62E+01 0.000 0.5132 0.9000 0.9000   1.00  1  1  9.5E-01
  6 :   9.62E-01 1.32E+01 0.000 0.5044 0.9061 0.9000   1.00  1  1  9.1E-01
  7 :   1.10E+00 8.32E+00 0.000 0.6297 0.9000 0.9484   1.00  1  1  8.7E-01
  8 :   1.16E+00 5.01E+00 0.000 0.6024 0.9000 0.6707   1.00  1  1  8.0E-01
  9 :   1.25E+00 1.93E+00 0.000 0.3861 0.9247 0.9000   1.00  1  1  6.0E-01
 10 :   1.32E+00 1.14E+00 0.000 0.5890 0.9000 0.9290   1.00  1  1  4.5E-01
 11 :   1.36E+00 7.26E-01 0.016 0.6373 0.9000 0.9450   1.00  1  1  3.4E-01
 12 :   1.38E+00 3.01E-01 0.000 0.4153 0.9255 0.9000   1.00  1  1  1.8E-01
 13 :   1.39E+00 1.76E-01 0.000 0.5839 0.9000 0.9327   1.00  1  1  1.1E-01
 14 :   1.40E+00 1.09E-01 0.051 0.6170 0.9000 0.9442   1.00  1  1  7.1E-02
 15 :   1.40E+00 5.69E-02 0.000 0.5238 0.9088 0.9000   1.00  1  1  3.8E-02
 16 :   1.40E+00 2.57E-02 0.000 0.4508 0.9000 0.9302   1.00  1  1  1.8E-02
 17 :   1.40E+00 1.04E-02 0.000 0.4070 0.9003 0.9000   1.00  1  2  7.2E-03
 18 :   1.40E+00 4.24E-03 0.000 0.4061 0.9000 0.9194   1.00  2  2  3.0E-03
 19 :   1.40E+00 1.40E-03 0.000 0.3292 0.9000 0.9096   1.00  2  2  9.7E-04
 20 :   1.40E+00 2.78E-04 0.000 0.1994 0.9000 0.9038   1.00  2  2  1.9E-04
 21 :   1.40E+00 3.73E-05 0.000 0.1339 0.9073 0.9000   1.00  2  2  2.6E-05
 22 :   1.40E+00 1.27E-06 0.000 0.0341 0.9904 0.9900   1.00  2  2  8.9E-07
 23 :   1.40E+00 2.67E-08 0.000 0.0211 0.9902 0.9900   1.00  2  2  1.9E-08
 24 :   1.40E+00 9.49E-09 0.000 0.3550 0.9344 0.9000   1.00  2  2  6.6E-09

iter seconds digits       c*x               b*y
 24      1.6   8.2  1.4016181622e+00  1.4016181527e+00
|Ax-b| =   4.3e-14, [Ay-c]_+ =   0.0E+00, |x|=  2.1e+00, |y|=  1.5e+03

Detailed timing (sec)
   Pre          IPM          Post
2.000E-02    1.590E+00    1.000E-02    
Max-norms: ||b||=1, ||c|| = 50,
Cholesky |add|=0, |skip| = 0, ||L.L|| = 1.99236.
------------------------------------------------------------
Status: Solved
Optimal value (cvx_optval): +1.40162
</pre>
<a id="plots"></a>
<div id="plotoutput">
<img src="l1_trend_filter_snp500__01.png" alt=""> 
</div>
</div>
</body>
</html>